Relaxation in a glassy binary mixture: Comparison of the mode-coupling theory to a 

Brownian dynamics simulation 



in 
o 
o 

< 



o 



Elijah Flenner and Grzegorz Szamel 
Department of Chemistry, Colorado State University, Fort Collms, CO 80525 

(Dated: February 2, 2008) 

We solved the mode-coupUng equations for the Kob- Andersen binary mixture using the structure 
factors calculated from Brownian dynamics simulations of the same system. We found, as was pre- 
viously observed, that the mode-coupling temperature, Tc, inferred from simulations is about two 
times greater than that predicted by the theory. However, we find that many time dependent quan- 
tities agree reasonably well with the predictions of the mode-coupling theory if they are compared 
at the same reduced temperature e= {T — Tc) /Tc, and if e is not too small. Specifically, the simula- 
tion results for the incoherent intermediate scattering function, the mean square displacement, the 
relaxation time and the self-diffusion coefficient agree reasonably well with the predictions of the 
mode-coupling theory. We find that there are substantial differences for the non-Gaussian parame- 
ter. At small reduced temperatures the probabilities of the logarithm of single particle displacements 
demonstrate that there is hopping-like motion present in the simulations, and this motion is not 
predicted by the mode-coupling theory. The wave vector dependent relaxation time is shown to 
be quahtatively different than the predictions of the mode-coupling theory for temperatures where 
hopping-like motion is present. 

PACS numbers: 64.70.Pf, 61.43.Fs, 61.20.Lc 
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I. INTRODUCTION 

The mode-coupling theory P, has been used exten- 
sively to describe the slow relaxation observed in super- 
cooled liquids close to the glass transition. The qualita- 
tive predictions of the theory have been compared to ex- 
periments 0,01 and simulations 0,|EIEI3- However, un- 
til recently, quantitative comparisons which include wave 
vector dependence and use exact (i.e. obtained from sim- 
ulations) static information have been rare H, 0, 0| . 

Nauroth and Kob Q compared the non-ergodicity pa- 
rameters predicted by the mode-coupling theory to New- 
tonian dynamics simulations of a binary Lennard-Jones 
mixture. They used the structure factors determined 
from the simulations as the input to the mode-coupling 
theory, and found that the transition temperature pre- 
dicted by the theory was approximately twice higher 
than the transition temperature inferred from simula- 
tions. However, the non-ergodicity parameters predicted 
by the theory agreed reasonably well with the ones deter- 
mined from simulations. In a later work, Kob, Nauroth, 
and Sciortino 9] quantitatively compared the shape of 
the intermediate scattering functions predicted by the 
mode-coupling theory to results of Newtonian dynamics 
simulations. Since the transition temperature predicted 
by the theory is greater than that inferred from simula- 
tions, the intermediate scattering functions at the same 
temperature (and close to the transition temperature) 
are vastly different. Therefore, the authors compared 
scattering functions predicted by the theory at a higher 
temperature to scattering functions obtained from sim- 
ulations at a lower temperature. The two correspond- 
ing temperatures were determined by the requirement 
that the relaxation time be the same for one wave vector 
around the first peak in partial structure factor for the 



larger particles. Once the two corresponding tempera- 
tures where fixed, the authors showed that for a few wave 
vectors the shape of the incoherent intermediate scat- 
tering function was well described by the mode-coupling 
theory. Also, Foffi et al. [l3 compared simulation re- 
sults to the mode-coupling theory for mixtures of hard 
spheres. They found that, if the time scale is rescaled, 
the mode-coupling theory accurately predicts the shape 
of the intermediate scattering functions in the a relax- 
ation region. 

Recently, the mode-coupling theory has been com- 
pared to molecular and Brownian dynamics simulations 
of polydisperse spheres with a strong repulsive core by 
Voigtmann, Puertes, and Fuchs They used the 

Percus-Yevick theory to determine the structure factors 
for the mode-coupling theory calculations. For the com- 
parison of the theory with simulation, they adjusted 
the packing fraction and allowed a "dynamical" length 
scale to vary slightly. They found that after these 
adjustments the intermediate scattering functions and 
the mean square displacement predicted by the mode- 
coupling theory agreed well with the simulation results. 

The goal of the work presented here is to compare the 
predictions of the mode coupling theory to the results of 
Brownian dynamics simulations using the smallest pos- 
sible number of adjustable parameters. In other words, 
we would like to test the predictive power of the mode- 
coupling theory rather than demonstrate that by making 
a number of parameters adjustable, we can reproduce the 
simulation results very accurately. In view of the differ- 
ence between the mode-coupling transition temperature 
predicted by the theory and inferred from simulations, a 
minimalistic approach is to compare theoretical predic- 
tions and results of the simulations at the same reduced 
temperature e = (T — Tc)/Tc. Briefly, the result of our 
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comparison is that many time dependent quantities pre- 
dicted by the theory agree reasonably well with results of 
the simulations if e is not too small. However, there is sig- 
nificant disagreement between the theoretical predictions 
and the simulation results for the non-Gaussian parame- 
ter for all but the highest reduced temperatures. Finally, 
at low reduced temperatures there is a hopping-like mo- 
tion present in the simulations which is not predicted by 
the mode-coupling theory. 

One of the reasons for testing the mode coupling the- 
ory against Brownian dynamics simulations is that the 
approximations involved in this theory applied to Brow- 
nian systems are somewhat less severe. The first ap- 
proximation of the theory, projecting stress fluctuations 
onto the subspace of the density products T], is exact 
for Brownian systems. However, the main and the most 
drastic approximation of the mode coupling theory, the 
self-consistent factorization approximation, has to still 
be used. The present study shows explicitly that it is 
the factorization approximation that is responsible for 
the failure of the mode coupling theory for low reduced 
temperatures. 

The paper is organized as follows. In section ^ 
we briefly describe the Brownian dynamics simulations. 
In section IIIII we describe the mode-coupling equations 
which are appropriate for a binary mixture evolving 
with Brownian dynamics, and briefly present the mode- 
coupling calculation (the method is described in detail 
in the appendix). In section Hvl we describe the method 
used to find the transition temperature from the mode- 
coupling theory, and briefly discuss the non-ergodicity 
parameters. We compare the intermediate scattering 
functions in section^and the mean square displacement 
in section IVII We compare the non-Gaussian parame- 
ter in section Ivnl and the probability distribution of the 
logarithm of single particle displacements in section lVIIII 
We discuss the results in section Hxl 
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with Vf being the gradient operator with respect to rf . 
In Eq. 1^ the random noise ffi satisfies the fluctuation- 
dissipation theorem. 



{Tj,{t)rf,{t'))^2DoS{t-t')S,,l. 
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In Eq. 0, the diffusion coefficient Do = fcsT/^o where 
ks is Boltzmann's constant and 1 is the unit tensor. 
Since the equation of motion allows for diffusive motion 
of the center of mass, all the results will be presented rela- 
tive to the center of mass (i.e. momentary positions of all 
the particles are always relative to the momentary posi- 
tion of the center of mass ) • The results are presented 
in terms of the reduced units with <jaa, ^aa, ^aaI^b^ 
and (j\a(x,I ^AA being the units of length, energy, tem- 
perature, and time, respectively. Since in these units the 
short-time self-diffusion coefficient is proportional to the 
temperature, in the comparisons with the mode-coupling 
theory the times are rescaled to t* — tDo/cr^A- 

The equations of motion, Eq. ^ were solved using a 
Heun algorithm with a small time step of 5 x 10^^. To 
save disk space, not all the generated configurations were 
saved to disk. Thus, the short time dynamics are not 
available at the lower temperatures. We simulated the 
temperatures T = 0.44, 0.45, 0.47, 0.50, 0.55, 0.60, 0.80, 
0.90, 1.0, 1.5, 2.0, 3.0, and 5.0. We ran equilibration 
runs and 4-6 production runs. The equilibration runs 
were typically twice shorter than the production runs, 
and the latter were up to 6 x 10^ time steps long for the 
lowest temperatures studied. The results presented are 
averages over the production runs. 



III. MODE-COUPLING THEORY 



II. BROWNIAN DYNAMICS SIMULATIONS 

We simulated a system consisting of Na = 800 par- 
ticles of type A and Nb = 200 particles of type B that 
was first considered by Kob and Andersen IH| . The inter- 
action potential is Va^(r) = 4eQ^[(CT„/3/r)i2 - {ac^p/rf], 
where a,^ S {A,B}, tAA = 1-0, oaa = 1-0, tAB = 1-5, 
CAB — 0.8, esB — 0.5, and obb = 0.88. The simula- 
tions are performed with the interaction potential cut at 
2.5 (Ta/3, and the box length of the cubic simulation cell 
is 9.4 gaa- Periodic boundary conditions were used. 

We performed Brownian dynamics simulations. The 
equation of motion for the position of the ith particle of 
type a, f f , is 
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where the friction coefficient of an isolated particle — 
1.0 and is the force acting on the ith particle of type 



The mode-coupling theory leads to a set of integro- 
differential equations for the coherent intermediate scat- 
tering fimctions (i.e. dynamic partial structure factors). 



5aM9,0 = (P?WP-,-(0) 



where 
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Note that the sum in Eq. Q is taken over particles 
of type a. The structure factors depend only on the 
magnitude of the wave vector \q \ = q. The time evolution 
of p2 (i) for a system of interacting Brownian particles is 



(1) governed by the adjoint Smoluchowski operator |l; 
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where f3 = l/fc^T and Dq is the short time diffusion 
coefficient which is the same for both types of particles. 
We set Dq = 1.0 for these calculations. 

The mode-coupling equations governing the time evo- 
lution of the coherent intermediate scattering functions 
for Brownian mixtures have been derived by Nagele et 
al. El: 



where the memory function 
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di 



S{q,t) = -q^DoS-\q)Siq,t) 



j duMiq,t-u) — Siq,u) 



where 



S(g,t) 



fSAA{q,t) SAB{q,t) 
\SBA{q,t) SBB{q,t) 
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is the matrix of coherent intermediate scattering func- 
tions, and M is the matrix of memory functions, 
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{\q-k\,t)Su'{k,t) (9) 



where the vortex 



Valmiq,k) = SamCalik) 



q- (g - fc) 



SclCa,m{\q~k\). (10) 



In Eq. IjlOl) the matrix C{q) is defined through the 
Ornstcin-Zernike matrix equation 



S-i(<z) = l-C(g), 



(11) 



where 1 is the unit tensor. The mode-coupling the- 
ory equations allows one to calculate the time evolu- 
tion of Sajsiq, t) with only the time independent quantity 
^apiq) = Sa/3iq,0) as an input. 

The incoherent intermediate scattering functions 



(12) 



are calculated using as input the coherent intermedi- 
ate scattering function and the partial structure factors. 
Nagele et al. derived the equations governing the 

time evolution of F^{q,t) for Brownian mixtures, 

^F^(g,i) = -q^D^F^^{q,t) 

duAf^((Z,t-u)|^i^^(q,u) (13) 



dk 



q ■ k 



Y,Co.5{k)Sss'{k,t)C^S'{k) (14) 



55' 



For short times the integrals involving the memory 
function in Eqs. lO and (|13(l are approximately zero, 
therefore 

S(g,i)«exp[-q2i^oS-i(<z)t] S(<z), (15) 



F'Jq,t)^e^j>[-q^Dot] . 



(16) 



In Eq. (O we used f ^(g,0) ^ 1.0. 

According to Eqs. (|15I16|) . at short times the parti- 
cles undergo diffusive motion with a diffusion coefficient 
Dq. The effect of the memory function is to provide a 
feedback mechanism which produces a "caging" of the 
particles. For temperatures below the transition temper- 
ature Tc, there is structural arrest and the particles do 
not escape their cage. This results in a non-zero value 
of S{q,t) and F^{q,t) as t — > oo. For temperatures close 
to but above Tc, there is a plateau region in the log- log 
plot of the mean square displacement and the log-linear 
plot of the intermediate scattering functions. At long 
times the motion of the particles is again diffusive with 
a temperature dependent diffusion coefficient D <C -Do- 

We calculate the input to the mode-coupling equations, 
i.e. the partial static structure factors Sap{q), directly 
from the Brownian dynamics simulations. Because of 
the finite size of the simulation box, the magnitude of the 
smallest wave vector calculated is Itt/L, where L is the 
length of the simulation box. We extrapolated Sap {q) to 
zero by fitting the first few wave vectors to a polynomial. 
The method used to calculate the integrals of the memory 
functions require that the partial structure factors are 
known at equally spaced wave vectors (see the appendix 
for details on the numerical procedure implemented in 
this work). The structure factor for these wave vectors 
were determined by fitting the partial structure factors 
determined from the simulations to a cubic spline. 

To solve the mode-coupling equations, we used 300 
equally spaced wave vectors from q = to g = 40 with 
the first wave vector go — 0.2/3. We performed a few 
calculations with larger cutoffs for the integral and/or a 
finer grid of wave vectors. The difference in the values 
of the calculated intermediate scattering functions was at 
most 5% and less in most cases. This difference results in 
a less than 5% difference in the self-diffusion coefficient 
and the incoherent intermediate scattering function's re- 
laxation time. 

The partial structure factors for temperatures which 
were not directly simulated were calculated by a 
quadratic polynomial interpolation between the points 
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at three adjacent temperatures. For most interpolated 
temperatures, it is possible to use two different sets of 
temperatures to determine the structure factor at the 
interpolated temperature. Using different temperature 
ranges changed the value of S{q,t) by as much as 1%, 
but by less in most cases. At a few temperatures we also 
used a linear interpolation between two adjacent tem- 
peratures or a cubic interpolation using the four closest 
temperatures. The difference in the results were less than 
2% using these different interpolation schemes. We con- 
clude that the results depends little on the interpolation 
scheme. 

First, we solved the mode-coupling equations for the 
coherent intermediate scattering functions until all the 
scattering functions decayed to zero or to a non-zero con- 
stant. Then the coherent scattering functions and the 
structure factor were used as the inputs for the calcu- 
lation of the incoherent scattering functions, which were 
solved for the same times as for the coherent intermediate 
scattering functions. 

Since the mode-coupling equations have to be solved 
for many decades in time, specialized techniques have 
been developed. We describe the method used to solve 
the mode-coupling equations in the appendix. The algo- 
rithm was first described by Fuchs et al. , and then, in 
more detail, by Miyazaki, Bagchi and Yethiraj It is 
an iterative technique which calculates the coherent and 
the incoherent scattering functions from t = ilSi where 
j G {A^ + 1, . . . , 2A^}, assuming that the scattering func- 
tions are known for t — jAt where j G {!,..., N}. Then 
the time step At' — 2At is doubled, and the values of the 
scattering functions for t = iAt' where i e {1, . . . , TV} are 
determined from the values of the scattering functions for 
t = jAt where j £ {1, . . . , 2A^}. We begin the calcula- 
tion for an initial time step of 10~*. For the first set of 
times, the scattering functions are not known. We used 
the approximation given by Eq. (|15|l and to supply 
the values of the coherent and the incoherent scattering 
functions for the initial N times. To check this approxi- 
mation, we also used a more time consuming procedure 
to solve the mode-coupling equations for the initial N 
times. In this procedure, the integrals of the memory 
function are included for t > 10~^. The difference in all 
calculated quantities was less than 0.5%. 



IV. TRANSITION TEMPERATURE: 
MODE-COUPLING THEORY VS. SIMULATIONS 



The mode-coupling theory predicts an ergodicity 
breaking transition when Sai3{q,t oo) ^ at some 
critical temperature Tc. To find the transition tempera- 
ture we calculate G{q) — limt^oo S(g,i) as a function of 
temperature. We followed the procedure used in previous 
calculations of the non-ergodicity parameter 0, |3, 0| to 
calculate G{q)- First, we took the Laplace transform of 
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FIG. 1: (a) Non-ergodicity parameter at Tc for the A parti- 
cles. Solid line: coherent correlator, A-A particles. Dashed 
line: incoherent correlator, A particles, (b) Non-ergodicity 
parameter at Tc'. coherent correlator, A-B particles, (c) 
Dashed line: A-B partial structure factor, SAB{q)- Solid line: 
GAsiq) = SAB{q;t oo) at Tc. (d) Non-ergodicity parame- 
ter at Tc for the B particles. Solid line: coherent correlator, 
B-B particles. Dashed line: incoherent correlator, B particles. 
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Eq.dto get 

[z + zM{q,z) + q^DoS-\q)] S(q, : 



S(g)+M(q,z)S(g) 
(17) 



where f{z) denotes the Laplace transform of f{t), f{z) = 
e~^*' fit). To find the long time limit of S{q,t) we 
utilized the relationship limt^oo f{t) = lini^^o ^fi^) and 
derived an equation for G{q). This equation was then 
solved using the following iterative procedure, 



q^Dol 



S(g)M«(q) 



S(q)M«(g)S(q), 
(18) 

where 'M.^{q) is calculated from G*(g). The non- 
ergodicity parameters are defined as fapil) = 
Gapiq)/Sap{q)- Using the same method, wc found the 
non-ergodicity parameters for the incoherent intermedi- 
ate scattering functions, F^{q,t oo) = fai^)- The 
iterative procedure for the non-ergodicity parameter for 
the self correlation functions is 



,(1+1) 



1 - q'Do ' 



(19) 



where Ma^'' depends on and the solution to Eq. 

The input to Eq. (|18|l are the partial structure fac- 
tors for a temperature T. The iterative procedure was 
followed until it was found that either G{q) is zero for 
all q or G{q) is nonzero and does not change anymore 
for all q. If G{q) is zero, then T > otherwise T < T^. 
By trying different structure factors for different temper- 
atures the transition temperature can be determined to 
arbitrary precision. This method is preferred to finding 
the transition temperature by calculating the full time 
dependence of S{q,t), since calculating G{q) is around 
20 times faster and reduces the calculation of the non- 
ergodicity parameter to hours instead of days. Once Tc 
was found, we calculated the non-ergodicity parameter 
by solving the mode-coupling equations for S{q,t) and 
F^(q,t) at Tc- The non-ergodicity parameters found us- 
ing the two methods agree. 

The non-ergodicity parameters at Tc are plotted in fig- 
ure ^ for the AA, AB, and BB correlators (solid lines), 
and for the incoherent intermediate scattering functions 
for the A and B particles (dashed lines). The results 
are similar to what was obtained by Nauroth and Kob 
0. All the non-ergodicity parameters are nonzero for 
q = 0, but approach zero for large q. The small fea- 
tures for q < 5.0 for the AA non-ergodicity parameters 
are numerical and do not represent additional features 
in the non-ergodicity parameter. The division of GAsiq) 
by SAsiq) for the AB non-ergodicity parameter causes 
numerical problems when SAsiq) ~ 0. This is seen as 
large spikes in f%g. In the insert we show the input to 
Eq. H18|l . SAsiq), (dashed line) and the results of the 



calculation Gab{<i) (solid line) at T^- An extensive com- 
parison of the mode-coupling theory predictions to the 
simulation results for the non-ergodicity parameters has 
already been conducted ^^l, and we do not repeat it here. 

We determined a transition temperature T*^'^"^'^ — 
0.9515, which is 3% higher than the ergodicity breaking 
temperature T^h.eory,NK ^ o.922 determined by Nauroth 
and Kob Foffi. et al. ^3 showed that small differ- 
ences in the structure factor can result in large differ- 
ences in the critical packing fraction for a binary system 
of hard spheres. They found that the critical packing 
fraction was around 5% higher if the partial structure 
factors were determined from the results of Newtonian 
dynamics simulations instead of using the Percus-Yevick 
approximation, even though there was little difference in 
the partial structures factors. 

The transition temperature predicted by the theory 
is around a factor of two larger than the commonly ac- 
cepted mode-coupling temperature inferred from simula- 
tions nil] of the same system, Tf" = 0.435. It should 
be emphasized that, in contrast to the transition pre- 
dicted by the mode-coupling theory, only a crossover in 
the dynamics is observed in simulations P, |3 ■ Fur- 
thermore, recently we argued that there is some arbitrari- 
ness regarding the mode-coupling temperature inferred 
from simulations [Toj . The mode-coupling temperature 
is usually obtained by fitting the simulation results for 
the characteristic decay time of the intermediate scatter- 
ing function and the diffusion coefficient to a power law 
a{T — Tc)'^ . The transition temperature obtained in this 
manner depends on the temperature range used in the fit. 
As argued in Ref. T^, around the commonly accepted 
value of the mode-coupling temperature, T|f '™ = 0.435, 
the relaxation mechanism changes from high temperature 
diffusive motion to low temperature hopping- like motion. 



V. INCOHERENT INTERMEDIATE 
SCATTERING FUNCTIONS 

We calculate the time dependence of the coherent and 
the incoherent intermediate scattering functions from the 
mode-coupling equations with the partial structure fac- 
tors determined from the simulations as input. In this 
section we compare predictions of the mode-coupling the- 
ory with results of Brownian dynamics simulations at the 
same reduced temperature e = [T — Tc) /T^- In tableQlwe 
list the reduced temperatures e, the corresponding tem- 
peratures for the Brownian dynamics simulations, and 
the temperatures which were used in the mode-coupling 
theory calculations. 

Incoherent intermediate scattering functions for the A 
particles are shown in Fig. and \^ calculated using 
the mode-coupling theory and obtained from the Brow- 
nian dynamics simulations, respectively. Both figures 
show the incoherent intermediate scattering function at 
q = 7.25, which is around the first peak of Saa{<i)- The 
scattering function for the mode-coupling theory was cal- 
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TABLE I: The reduced temperatures and their corresponding 
temperature in the Brownian dynamics simulations and the 
mode-couphng theory calculations. 



e 


BD Temperature 


MCT Temperature 


0.0000 


0.435 


0.9515 


0.0115 


0.44 


0.9624 


0.0345 


0.45 


0.9843 


0.0805 


0.47 


1.0281 


0.1495 


0.50 


1.0937 


0.2644 


0.55 


1.2030 


0.3793 


0.60 


1.3124 


0.8391 


0.80 


1.7499 


1.0690 


0.90 


1.9686 


1.2989 


1.00 


2.1873 


2.4483 


1.50 


3.2810 


3.5977 


2.00 


4.3747 


5.8966 


3.00 


6.5621" 


10.494 


5.00 


10.937" 



"mode-coupling equations were not solved for these temperatures. 



culated through a linear interpolation of nearby scatter- 
ing functions which were calculated directly from the 
mode-coupling equations. The dotted line in Fig. 12^- 
is the mode-couphng results at T^*'»'=°'-!' = 0.9515. The 
dashed line in each figure shows the incoherent scat- 
tering function for non- interacting particles, Eq. (|16|l . 
The mode-coupling theory correctly predicts the two step 
decay of the intermediate scattering functions, and the 
plateau observed in the log-linear plot of the scattering 
functions. 

In Fig. we compare the incoherent intermediate 
scattering functions for e = 3.5977, 0.8391 0.0805, and 
0.0115 for the A particles a,t q — 7.25, and in Fig.|3jD we 
show the comparison for the same reduced temperatures 
for the B particles at g = 5.75. At the higher reduced 
temperatures there is very good agreement between the 
mode-coupling calculations and the Brownian dynamics 
simulations. For 0.0345 < e < 0.8391 the characteristic 
decay time of the scattering functions calculated from the 
mode-coupling theory is less than that of the Brownian 
dynamics simulation. For reduced temperatures equal 
to and below 0.0345, the mode-coupling theory predicts 
a longer decay time for the self intermediate scattering 
functions for this value of q. However, the shape of the 
incoherent intermediate scattering functions arc similar 
in the a relaxation region. 

Kob, Nauroth, and Sciortino |s|| compared the self 
intermediate scattering functions obtained from Newto- 
nian dynamics simulations and predicted by the mode- 
coupling theory. The input temperature in the mode- 
coupling theory was adjusted so that the relaxation time 
of SAA{q,t)/ SAA{q) was correctly reproduced for one 
wave vector around the first peak of SAA{q)- The proce- 
dure followed by Kob et al. requires that the character- 
istic decay time is the same in the simulations and the 
mode-coupling theory calculations. Kob et al. observed 
that the shape of the scattering functions and its wave 





toy 

FIG. 2: Incoherent intermediate scattering function for the A 
particles for q — 7.25. (a) Predicted by the mode-coupling 
theory, (b) Calculated from the Brownian dynamics sim- 
ulations. The scattering functions are shown for the same 
reduced temperatures e = (T — Tc)/Tc. The reduced tem- 
peratures are 3.5977, 2.4483, 1.2989, 1.0690, 0.8391, 0.3793, 
0.2644, 0.1494, 0.0805, 0.0345, 0.0115 listed from left to right. 
The dashed line in both figures correspond to the limit of non- 
interacting particles. The dotted line in figure (a) is the inco- 
herent intermediate scattering function calculated at 



vector dependence is accurately described by the theory 
for reduced temperatures above 0.071. We discuss the 
wave vector dependence of the relaxation time in section 

The mode-coupling theory predicts power law diver- 
gence of the characteristic decay time of the intermedi- 
ate scattering functions. Specifically, we define the a 
relaxation time as the time when the incoherent inter- 
mediate scattering function decays to of its initial 
value, F^{q,Ta) — e^^. In Fig. Q we show the a re- 
laxation time for the A and B particles as a function of 
reduced temperature. We fit the a relaxation time to 
the function a[{T - Tc)/Tc]^^. We fit the Brownian dy- 
namics results to reduced temperatures from 0.8391 to 
0.1495, which corresponds to the same reduced temper- 
atures fit in an earlier work We used this range 
of temperatures since this is the temperature range in 
which a power law fits the a relaxation time well for a 
transition temperature of 0.435. For the mode-coupling 
theory calculations, we fit the power law to reduced tem- 
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FIG. 3: The incoherent intermediate scattering function for 
the A and B particles predicted by the mode-coupling the- 
ory (dashed lines) and calculated from the Brownian dynam- 
ics simulations (solid lines) for e — 3.598, 0.839, 0.0805 and 
0.0115 listed from left to right. 



peratures equal to and below 0.0345. The exponents in 
the power law fits are given in the figure. 

The exponents from the mode-coupling calculations 
are close to but slightly larger than the exponents found 
from the simulations. Note that the exponent obtained 
from solving the mode-coupling equations for the single 
component hard sphere system is the same as the ex- 
ponents predicted by the mode coupling coupling theory 
for the binary Lennard- Jones system. The exponents de- 
termined from simulations are slightly different than the 
exponents reported in an earlier work p^ . since we are 
fitting TqDq here and we fit Tq in the other work. We 
find, as other authors have also observed that the 
exponents determined from the simulations are different 
for the A and the B particles. However, this difference is 
within the 3% uncertainty in the exponents. Note that 
the mode coupling theory predicts the same exponents 
for the A and B particles, and the same exponents for 
the a relaxation time and the self diffusion coefficient. 

Furthermore, the predictions of the mode-coupling the- 
ory follow the asymptotic power law behavior for reduced 
temperatures smaller than e ~ 0.08. This is in contrast 
with the results of simulations, which can only be fit to 
power laws in a restricted range 0.1495 < e < 0.8391. 
On the other hand, it has been observed in experiments 
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FIG. 4: The a relaxation time calculated from the Brownian 
dynamics simulations (symbols) and predicted by the mode- 
coupling theory (solid and dashed lines) for the the A and 
B particles. The A particles are represented by the closed 
symbols and the solid line. The B particles are represented 
by the open symbols and the dashed line. The lines nearly 
overlap. The dotted lines are fits of the simulation data to 
the function a[{T — Tc) /Tc]~'^ ■ The exponents to the power 
law fits are given in the figure along with the exponents to 
power law fits to the predictions of the mode-coupling theory. 



that the power law behavior of Tq is valid for reduced 
temperatures up to 0.5 poj |. 



VI. MEAN SQUARE DISPLACEMENT 

To derive the mode-coupling theory predictions for the 
mean square displacement we use the small q expansion 
of the incoherent intermediate scattering function 23], 

F^t) = ^(-1)»^_^ (Srl-it)) , (20) 

where {Sr^J'it)) = (|rf (t) - rf (0)|2"). Inserting 
Eq. H20|l into the equation for the incoherent intermedi- 
ate scattering function, Eq. p3|l . and expanding M^(q, t) 
in a Taylor series 

M^(q, t) = Mlit) + q'Mlit) + q'Mtit) + ... (21) 

results in the equation of motion 

^^{5rl{t))^&D,^ j\uMl{t~u)^{5rl{u)) (22) 



where 



Ml{t) = \mimq.t) 



dfc k^F^{k,t) 



J2Cas{k)Sss'{k,t)C^S'{k). (23) 
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FIG. 5: Mean square displacement for the A particles, (a) 
Predicted by the mode-coupling theory, (b) Calculated from 
the Brownian dynamics simulations. The solid lines corre- 
spond to the same reduced temperatures e = {T — Tc)/Tc 
in the mode-coupling theory calculations and the Brownian 
dynamics simulations. The reduced temperatures are 3.5977, 
2.4483, 1.2989, 1.0690, 0.8391, 0.3793, 0.2644, 0.1494, 0.0805, 
0.0345, 0.0115 listed from left to right. The dashed line corre- 
sponds to the limit of non-interacting particles and the dotted 
line in figure (a) is the mean square displacement calculated 
at y*'"^"'"!' 



This equation was solved using a variation of the proce- 
dure described in the appendix. 

We present the mode-coupling theory predictions and 
simulations results in Figs. and[S|3, respectively, for 
the A particles. The results for the B particles are sim- 
ilar. The solid lines correspond to the same reduced 
temperatures for the mode-coupHng calculations and the 
Brownian dynamics simulations. The dotted line in 
Fig. 01 is the mean square displacement at Tl^'^°'^y . The 
dashed line in both figures corresponds to motion in the 
limit of non- interacting particles, i.e. a purely diffusive 
motion with a diffusion coefficient Dq. At all tempera- 
tures the short time motion is diffusive with a diffusion 
coefficient Dq. For the mode-coupling theory calcula- 
tions, the short time diffusion coefficient is one, but the 
short time diffusion coefhcient is temperature dependent 
in the Brownian dynamics simulations. Note, however, 
that the time axes are scaled to compensate for this dif- 
ference. 




FIG. 6: The mean square displacement for the A and B par- 
ticles predicted by the mode-coupling theory (dashed lines) 
and calculated from the Brownian dynamics simulations (solid 
lines). The reduced temperatures are e = 3.5977, 0.839, 
0.0805 and 0.0115 listed from left to right. 



The mode-coupling theory correctly predicts the ex- 
istence of the plateau region in the log-log plot of the 
mean square displacement for temperature close to the 
transition temperature. The plateau represents a local- 
ization of the particles for several decades in time and is 
associated with the cage effect. According to the mode- 
coupling theory, above the transition temperature the 
long time motion is diffusive with a self diffusion coeffi- 
cient D > Q. Also, at and below the mode-coupling tran- 
sition temperature, there is structural arrest and D — Q. 

In Fig. we show the mean square displacements for 
reduced temperatures e = 3.5977, 0.8391 0.0805, 0.0115 
calculated using the mode-coupling theory (dashed lines) 
and obtained from the Brownian dynamics simulations 
(solid lines) for the A and B particles. The A particles 
are shown in the upper figure and the B particles are 
shown in the lower figure. The mean square displace- 
ment agrees reasonably well with the predictions of the 
mode-coupling theory for reduced temperatures equal to 
and above 0.0805 for both the A and B particles at all 
times. However, for reduced temperatures below 0.0805, 
the mean square displacement is not accurately predicted 
by the mode-coupling theory. 

For all reduced temperatures smaller than e — 0.8391, 
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FIG. 7: Comparison of the plateau value of the mean square 
displacement versus temperature. The symbols are the sim- 
ulation results and the lines are the predictions of the mode- 
coupling theory. The A particles are represented by the closed 
symbols and the solid line, and the B particles are represented 
by the open symbols and the dashed line. 



there is an obvious plateau in the Brownian dynam- 
ics simulations and in the mode-coupling theory calcu- 
lations. We define the height of the plateau region as 
the inflection point in the logarithm of the mean square 
displacement versus the logarithm of time. At small re- 
duced temperatures the height of the plateau predicted 
by the theory agrees reasonably well with that obtained 
from simulations, Fig. In particular, at e = 0.0115, 
the plateau is predicted by the theory at a value of the 
mean square displacement of around 0.028 for the A par- 
ticles and it occurs in the Brownian dynamics simula- 
tions at around 0.029. For the B particles at e = 0.0115, 
the plateau is predicted by the theory to be around a 
value of 0.053 and it is around 0.043 in the simulations. 
However, the temperature dependence of the plateau ac- 
cording to the theory and in the simulations is very dif- 
ferent. The theory predicts that the plateau height as 
a function of reduced temperature is essentially constant 
until around e = 0.8391, and the plateau height increases 
slightly with increasing temperature. The plateau height 
calculated from the simulations increases with temper- 
ature faster than predicted by the theory. For reduced 
temperatures above 0.3793, it is difficult to calculate the 
infiection point for the Brownian dynamics simulations 
accurately. Note that in the temperature range in which 
the theory gives reasonably accurate predictions for the 
incoherent scattering function and the mean square dis- 
placement, the plateau height resulting from the theory is 
quite a bit smaller than that obtained from simulations. 
For example, for e = 0.3793 the mean square displace- 
ment at the inflection point predicted by the theory is 
around 0.024 and 0.053 for the A and B particles, re- 
spectively, whereas in simulations it occurs around 0.040 
and 0.12 for the A and B particles, respectively. 

We obtain the self diffusion coefficient D from the slope 
of the mean square displacement at long times, i.e. from 
D = limi_oo /(6i)- Within the mode-couphng 
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FIG. 8: The diffusion coefficient determined from the Brow- 
nian dynamics simulation (symbols) and the mode-coupling 
theory (solid and dashed lines). The closed symbols and the 
solid line are the results for the A particles. The open sym- 
bols and the dashed line are the results for the B particles. 
The dotted line is a fit of the simulation data to the function 
a\{T — T^ jT^ . The exponents to the power law fits are given 
in the figure along with the exponents to power law fits to the 
predictions of the mode-coupling theory. 



theory, we can also calculate D from the equation 
D _ 1 



(24) 



Both procedures agree to within one percent. The 
diffusion coefficients predicted by the theory and ob- 
tained from Brownian dynamics simulations are shown in 
Fig. IHlas a function of reduced temperature. The mode- 
coupling theory provides a good prediction for the diffu- 
sion coefficients for e > 0.0805, but the diffusion coeffi- 
cients predicted by the theory are slightly larger than the 
ones found from simulations. At lower reduced temper- 
atures the theory strongly underestimates the diffusion 
coefficients. In addition, the theory does not capture the 
increasing difference between the diffusion coefficients of 
the A and the B particles. 

We fit D/Do to power laws of the form a{[T-T^) jT^^; 
the exponents are given in Fig.|Sl For the mode-coupling 
theory fit we use reduced temperatures less than 0.0345 
whereas for the Brownian dynamics simulations we fit 
the range 0.1495 < e < 0.8391. The exponents for 
the Brownian dynamics simulations are slightly different 
than what has been reported in an earlier work |19l] , since 
we are fitting D/Dq here and D there. The exponents 
predicted by the mode-coupling theory are considerably 
larger than those obtained from the fits to simulations' 
results. Note that the exponents for the diffusion coef- 
ficient calculated from the theory is the same as the ex- 
ponents found for the a relaxation time. Moreover, the 
exponents found from the mode-coupling calculations are 
the same for the A and the B particles. 

The reduced temperatures in which it is possible to 
fit the mode-coupling results well with a power law is 
similar to what was found for the a relaxation time, i.e. 
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FIG. 9: Product of the diffusion coefficient and the a re- 
laxation time determined from the Brownian dynamics sim- 
ulation (symbols) and the mode-coupling theory (solid and 
dashed lines). The closed symbols and the solid line are the 
results for the A particles. The open symbols and the dashed 
line are the results for the B particles. Note that for clarity 
we omitted the error bars for the lowest temperature point 
(e — 0.0115) for the A particles; at the lowest temperature 
for the A particles Dtc, = 0.0585 ± 0.054 



that the power law provides a good fit up to a reduced 
temperature of around 0.08. 

We would like to point out that power laws fit the 
predictions of the mode-coupling theory reasonably well 
for the same reduced temperatures in which we fit power 
laws to the results of the Brownian dynamics simulations. 
If we fit the predictions of the mode-coupling theory to 
power laws using the range 0.1495 < e < 0.8391, the re- 
sulting exponents are different from the true exponents 
{i.e. from the exponents describing the the true asymp- 
totic power law behavior) but they differ by at most 
12% from the exponents obtained from the Brownian dy- 
namics simulations. This should be compared with the 
48% difference between the true mode-coupling theory 
exponent and the exponent for the B particles obtained 
from the Brownian dynamics simulations using the range 
0.1495 < e < 0.8391. 

Finally, in Fig. I^we compare the temperature depen- 
dence of the product of the diffusion coefficient and the 
a relaxation time predicted by the mode-coupling theory 
and obtained from simulations. In the high temperature 
regime one usually expects that the Stokes-Einstein rela- 
tion is valid and Dtq, is temperature- independent. The 
decoupling of the diffusion and structural relaxation has 
been identified as one of the signatures of increasingly 
heterogeneous dynamics |2ll | . The mode-coupling theory 
predicts an essentially temperature-independent Dra . In 
contrast, simulations show that the product of the diffu- 
sion coefficient and the a relaxation time starts increasing 
with temperature below approximately e = 1. Note that 
this value of the reduced temperature corresponds to a 
temperature that is close to the so-called onset temper- 
ature identified for the Kob- Andersen model by Brumer 
and Reichman (22i |. 



VII. NON-GAUSSIAN PARAMETER 

At short and long times, the motion of the particles is 
Fickian and the self part of the van Hove correlation func- 
tion is Gaussian. For intermediate times, the van Hove 
correlation function deviates from Gaussian. To exam- 
ine the non-Gaussian nature of the van Hove correlation 
function, we calculated the non-Gaussian parameter 



a2{t) 



3 (Srtit)) 
5 (Srlit))' 



(25) 



For the calculation of a2 (t) , we first calculated the mean 
square displacement (see section IVI|I and then we cal- 
culated {Sr'^{t)'). The equation of motion for (Sr^it)) 
is derived using the same method as the equation for 
{5r1^{t)). The resulting equation of motion is 



2QDo{Srl{t)) 

+ 10 j\uMl{t~u)^{Srl{u)) 



(26) 



where 
Mlit) 



107r2Ar„ 
dfc k'^ 



2d d'^ 

K{k,t) + —F-^{k,t) 



Y,Cas{k)Sss'{k,t)Co^S'{k). 



(27) 



SS' 



In principle, it is also possible to calculate the non- 
Gaussian parameter by fitting F^(q, t) in the small wave 
vector limit We found that this procedure was 

difficult to follow using the structure factors calculated 
from simulations, and small numerical uncertainties in 
the small q values of can result in large changes in the 
non-Gaussian parameter. For short times, the integrals 
involving memory functions in Eqs. H22|l and 1)26(1 are 
close to zero. Thus, for short times {Sr^{t)') « 6Dot and 
(Sr'^it)) « GOD^t"^, therefore 02(1) « 0. At long times, 
the motion is once again Fickian, the self part of the van 
Hove correlation function is Gaussian and 02 (i) = 0. In 
our calculations, a2{t) does not go to zero at long times, 
but rather to a value around 5 x 10~^. We believe that 
this result can be attributed to numerical errors; it does 
not imply that a2{t) is nonzero for t — > 00. 

We show the non-Gaussian parameters for the A par- 
ticles in Fig. El The upper panel shows the predictions 
of the mode-coupling theory and the lower panel shows 
the results of the Brownian dynamics simulations. Ac- 
cording to the mode-coupling calculations, there is one 
peak at high reduced temperatures, but for e < 1.2989 
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FIG. 10: The non-Gaussian parameter for the A particles, 
(a) Predicted by the mode-coupling theory, (b) Calculated 
using the Brownian dynamics simulations. The reduced tem- 
peratures e = (T - T^)/Tc are 3.5977, 2.4483, 1.2989, 1.0690, 
0.8391, 0.3793, 0.2644, 0.1494, 0.0805, 0.0345, 0.0115. In fig- 
ure (a), there is one peak for e = 3.5977, and a wider single 
peak for 2.4483 (dashed lines). For all the other reduced tem- 
peratures there are two peaks, and the larger peak position 
of the second peak corresponds to a lower reduced temper- 
ature. In figure (b), the larger peak heights correspond to 
lower reduced temperatures. 



there are two peaks. Two peaks have been observed be- 
fore in mode-coupUng calculations of a2 (t) for other sys- 
tems 23;] . The first peak is around the beginning of the 
plateau region of the mean square displacement and the 
initial decay of the scattering functions. The position of 
the first peak decreases for decreasing temperature, but 
is almost constant for e < 0.8391. The second peak is 
around the a relaxation time which corresponds to just 
after the plateau region of the mean square displacement. 
The position of the second peak increases with decreasing 
temperature and roughly follows the temperature depen- 
dence of the a relaxation time, Fig. 112b . 

There is only one peak in a2 {t) according to the Brow- 
nian dynamics simulations, Fig. llOb . The position of this 
peak is greater than the a relaxation time at higher tem- 
peratures, but it increases slower with decreasing tem- 
perature than the a relaxation time starting at around 
e = 0.8391. 

The mode coupling theory predicts one peak at all tem- 



FIG. 11: The non-Gaussian parameter for the B particles, 
(a) Predicted by the mode-coupling theory, (b) Calculated 
from the Brownian dynamics simulations. The reduced tem- 
peratures e = (T - T^)/Tc are 3.5977, 2.4483, 1.2989, 1.0690, 
0.8391, 0.3793, 0.2644, 0.1494, 0.0805, 0.0345, 0.0115 where 
the lower reduced temperatures correspond to larger peak 
heights. 



perature for the B particles. However, there is a promi- 
nent shoulder for e < 1.2989, the same reduced temper- 
atures in which there are two peaks in the non-Gaussian 
parameter for the A particles. The position of the shoul- 
der follows the same temperature dependence as the po- 
sition of the first peak in the non-Gaussian parameter 
for the A particles. The peak position in a2{t) predicted 
by the mode-coupling theory is less than the a relaxation 
time at all temperatures, and it increases with decreasing 
temperature at close to the same rate as the a relaxation 
time, Fig, nib. 

The non-Gaussian parameter for the B particles ob- 
tained from the Brownian dynamics simulation is shown 
in Fig.lllb. There is only one peak for all temperatures. 
The position of the peak is greater than the a relaxation 
time for higher temperatures, but the position increases 
slower with decreasing temperature than the a relaxation 
time and is much less than the a relaxation time at small 
reduced temperatures. Fig. WBa . 

To summarize, the time dependence of the non- 
Gaussian parameter calculated using the mode-coupling 
theory is significantly different from what is obtained 
from the Brownian dynamics simulations. More im- 
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FIG. 12: The peak position of the non-Gaussian parameter 
compared to the a relaxation time. The symbols are the sim- 
ulation results and the lines are the predictions of the mode- 
coupling theory. The open symbols and dashed lines are the 
a relaxation time. The closed symbols and the solid lines are 
the peak positions of the non-Gaussian parameter. 



portantly, the mode-coupling theory strongly underesti- 
mates the deviations from Gaussian {i.e. Fickian) dif- 
fusive motion: the heights of non-Gaussian parameters 
predicted by the theory are almost an order of magnitude 
smaller than those obtained from the Brownian dynam- 
ics simulations. We show in the next section that the 
simulations show even stronger non-Gaussian effects on 
somewhat longer time scales. 



VIII. PROBABILITY DISTRIBUTIONS OF 
THE LOGARITHM OF SINGLE PARTICLE 
DISPLACEMENTS 

In an earlier work we showed that as the tem- 
perature = 0.435 is approached, the motion of the 
particles changes from a high temperature diffusive-like 
behavior to a low temperature hopping-like motion. To 
this end, we investigated probability distributions of the 
logarithm of single particle displacements 0, ^| . The 
probability distribution of the logarithm of single par- 
ticle displacements at a time t, P{logiQ{Sr);t), can be 
obtained from the self van Hove correlation function, 
GsiSr.t) = {Si\ff{t)-ff{0)\-Sr)), by the following 




FIG. 13: The probability distribution of the logarithm of 
single particle displacements predicted by the mode-coupling 
theory (dotted hues) compared to the probability distribu- 
tions calculated from the Brownian dynamics simulations 
(solid lines) for e = 0.8391. (a) A particles, (b) B parti- 
cles. The times shown are 0.25, 1.0, 10, and 20 times the 
a relaxation time; note that we use the A-particles and the 
B-particles a relaxation time in (a) and (b), respectively. 



transformation P(logig((5r); t) = ln{10)A'KSr^GsiSr,t). 
Note that if the motion of a tagged particle is diffu- 
sive at all times with a diffusion coefficient D, then 
the self van Hove correlation function Gg {Sr, t) = 
(l/(47r£'t) 2 ) exp(— 5r^/4_Di) and it can be shown that 
the shape of the probabihty distribution P{\ogiQ{Sr);t) 
is time-independent. In particular, the peak height of 
P{\ogiQ{Sr);t) does not depend time and is equal to 
ln(10)v/54/7r e"3/2 _ 2.13. Thus, deviations from this 
height represent deviations from Gaussian behavior of 
G,((5r,t). 

We show a comparison of P{\ogif^{Sr);t) in Fig. 1131 for 
the A and B particles calculated from the simulations 
(solid lines) and from the mode-coupling theory (dotted 
lines) for the reduced temperature e = 0.8391 for several 
different times: 0.25, 1.0, 5.0, and 10.0 Tq. It should 
be noted that in Figs. 1131 and ITU we use the A and the 
B particles a relaxation times in panels (a) and (b), re- 
spectively; moreover, we use the a relaxation times ob- 
tained from simulations for the simulation results and the 
a relaxation times predicted by the mode-coupling the- 
ory for the theoretical results. At this temperature the 
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FIG. 14: The probability distribution of the logarithm of 
single particle displacements predicted by the mode-coupling 
theory (dotted lines) compared to the probability distribu- 
tions calculated from the Brownian dynamics simulations 
(solid lines) for e = 0.0115. (a) A particles, (b) B parti- 
cles. The times shown are 0.1, 1, 2.5, 5, and 10 times the 
a relaxation time; note that we use the A-particles and the 
B-particles a relaxation time in (a) and (b), respectively. 



mode-coupling theory describes the probability distribu- 
tions reasonably well. 

In Fig. E|we show a comparison of P (\ogyJy8r)\t^ for 
the A and B particles calculated from the simulations 
(solid lines) and from the mode-coupling theory (dotted 
lines) for several different times at a reduced tempera- 
ture of e = 0.0115 (recall that, as in Fig. 1131 we use the 
A and the B particles a relaxation times in panels (a) and 
(b), respectively; moreover we use the a relaxation times 
obtained from simulations for the simulation results and 
the a relaxation times predicted by the mode-coupling 
theory for the theoretical results). This is the lowest re- 
duced temperature in which we can directly compare the 
predictions of the mode-coupling theory to the simula- 
tion results. There is a dramatic difference in the shape 
of the curves over the time interval shown in the fig- 
ure. At intermediate times bimodal distributions are ob- 
tained from Brownian dynamics simulations whereas the 
mode-coupling theory predicts unimodal distributions at 
all times. The bimodal distributions suggests that a por- 
tion of the particles are undergoing hopping-like motion 
with a large distribution of hopping rates. This hopping- 




FIG. 15: The product D<jVa(<j) for e = 0.8391 and 0.0115. 
The symbols are the results of the Brownian dynamics simu- 
lations and the lines are the predictions of the mode-coupling 
theory. The A particles are represented by the closed symbols 
and the solid lines, and the B particles are represented by the 
open symbols and the dashed lines. 



like motion is not predicted by the mode-coupling theory. 

In an earlier work , we observed in the simulations 
that the time in which both peaks are about equal height 
is longer than the time indicated by the peak position of 
the non-Gaussian parameter a^it)- We defined a new 
parameter ^{t) — i ((5r^) (l/(5r^) — 1 whose peak po- 
sition occurs at the same time as when the two peaks 
in P(log]^Q((5r); are approximately the same height for 
temperatures in which there are two peaks. It was found 
that the peak position of 7(t) has the same temperature 
dependence as the a relaxation time. 

At temperatures in which we see evidence of the 
hopping-like motion when we examine the probability 
distribution P(logiQ((5r); i), the wave vector dependent 
a relaxation time is qualitatively different in the simula- 
tion than predicted by mode-coupling theory. 

In Fig. El we show Dq^Ta{q), where Ta{q) is the wave 
vector dependent a relaxation time defined as the time 
when F^{q, Ta{q)) — e~^. The lines are the predictions of 
the mode-coupling theory and the circles are the results 
of the Brownian dynamics simulations. For small wave 
vectors, the product Dq^Ta{q) is one. For the reduced 
temperature of e = 0.8391, Fig. [T3i, T„(g)~^ = Dq^ 
until q « 5. At larger wave vectors there is a crossover 
from the small wave vector relationship to the large wave 
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vector limit Ta{q)~^ = Dgq^ (this large wave vector limit 
follows from the fact that memory functions vanish in the 
large wave vector limit). 

In contrast, for lower reduced temperatures, there is a 
qualitative difference between the wave vector dependent 
a relaxation time predicted by the theory and calculated 
from the Brownian dynamics simulations. In Fief. ll5b we 
show the product Dq^Ta{q) for the reduced temperature 
e = 0.0115. At small wave vectors the simulation results 
approach the asymptotic behavior Ta{q)~^ = Dq^. At 
intermediate wave vectors, 1.5 < q < 7.5, the values of 
Dq^Ta{q) calculated from the simulation increase with a 
peak somewhere between 6.5 < q < 7.5. At large wave 
vectors, the a relaxation time reaches its limiting behav- 
ior Ta{q)~'^ = Doq^. On the other hand, the mode cou- 
pling theory predicts a monotonically decaying Dq^Ta{q) 
with essentially the same behavior at low and high re- 
duced temperatures. Fig 1 151 

We should point out that the increase of Dq^Ta{q) with 
wave vector for intermediate wave vectors was previously 
found in the Kob- Andersen system by Berthier [2g. Also, 
the increase of Dq^Ta{q) has been predicted within the 
theoretical approach of Schweizer and Saltzman |23| . Fi- 
nally, Berthier, Chandler and Garrahan found a similar 
behavior in one-dimensional facilitated kinetic Ising mod- 
els |2|| . Here we emphasize that this behavior correlates 
with strong deviations from Fickian diffusion visible in 
the probability distribution P(logiQ((5r); t). 



IX. CONCLUSIONS 

We have conducted an extensive comparison of the pre- 
dictions of the mode-coupling theory to Brownian dy- 
namics simulations. As has been previously observed, 
qualitatively, predictions of the mode-coupling theory 
agree well with simulations. Namely, the mode-coupling 
theory accurately predicts the two step relaxation ob- 
served in the scattering functions, the existence of the 
plateau region observed in the mean square displacement, 
and the power law behavior of the self diffusion coefficient 
and the a relaxation time. 

The mode-coupling theory overestimates the feedback 
mechanism in the memory functions, resulting in a tran- 
sition temperature T*'^^"^'^ = 0.9515 much greater than 
the temperature T^*™ = 0.435 inferred from simulations. 
The transition temperature found in simulations is de- 
termined by fitting the diffusion coefficient and the a 
relaxation time to power laws. It should be noted that 
this temperature is sensitive to the range of tempera- 
tures used in the power law fits While the transi- 
tion temperatures are vastly different, we find that the 
mode-coupling theory gives good quantitative results of 
many time dependent quantities if they are compared at 
the same reduced temperature e — {T — Tc)/Tc. 

The self intermediate scattering functions and the 
mean square displacement resulting from simulations are 
well described by the mode-coupling theory for reduced 



temperatures greater than 0.08. For temperatures close 
to Tc the mode-coupling theory predicts divergence of 
the self intermediate scattering function's relaxation time 
and vanishing of the self diffusion coefficient. There is no 
divergence of the relaxation time or vanishing of the self 
diffusion coefficient in the Brownian dynamics simula- 
tions. 

The non-Gaussian parameter calculated from the sim- 
ulations are quite different than the non-Gaussian param- 
eter predicted by the mode-coupling theory. The mode- 
coupling theory predicts two peaks in a2{t) for reduced 
temperatures e < 1.2989 for the A particles, and a shoul- 
der at short times and a peak at longer times for the B 
particles. There is only one peak in the non-Gaussian 
parameter calculated from the simulation for all tem- 
peratures for both the A and B particles. Furthermore, 
the position of the second peak for the A particles and 
the only peak for the B particles predicted by the mode- 
coupling theory has a different temperature dependence 
than the position of the peak observed in the Brownian 
dynamics simulations. The mode-coupling theory pre- 
dicts that the position of the second peak roughly fol- 
lows the a relaxation time close to Tc, while the position 
of the peak for the Brownian dynamics simulations in- 
creases slower with decreasing temperature than the a 
relaxation time. Finally, the theory underestimates the 
height of the non-Gaussian parameter peak by almost an 
order of magnitude. 

We calculated the probability of the logarithm of single 
particle displacements P{logiQ{5r);t) which is sensitive 
to hopping-like motion. At high reduced temperatures, 
there is no hopping-like motion evident in the Brow- 
nian dynamics simulations and P{\ogiQ (dr);t) is accu- 
rately described by the mode-coupling theory. For low re- 
duced temperatures, there is little agreement between the 
mode-coupling theory and the simulations. At the lowest 
reduced temperature studied in this work, e — 0.0115, 
the hopping-like motion is very evident in the Brown- 
ian dynamics simulations, but is not predicted by the 
mode-coupling theory. We believe that this hopping-like 
motion is responsible for the absence of the divergence 
of the relaxation time or vanishing of the self diffusion 
coefficient in the simulations. 

For low reduced temperatures, the mode-coupling the- 
ory does not predict the proper wave vector dependence 
of the a relaxation time. The theory predicts that at all 
temperatures the product Dq^T{q) is one for small wave 
vectors and it decreases monotonically with increasing 
wave vector. In contrast, at low temperatures Brown- 
ian dynamics simulations show a peak in the product 
i^gV(g). 

To summarize, we found that close but not too close 
to the transition temperature the mode-coupling theory 
does predict most time-dependent quantities reasonably 
well. At very small reduced temperatures there is a 
hopping-like motion present in the simulations which is 
not accounted for by the standard version of the mode- 
coupling theory. Signatures of the hopping-like motion 
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include the two-peak structure of the probabihty of the 
logarithm of particle displacements and non-trivial wave 
vector dependence of the a relaxation time. 
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APPENDIX: NUMERICAL ROUTINES 

In this section we will outline the numerical routines 
that were used to solve many of the equations given in 
this paper. Since the time dependent quantities were 
solved over many decades in time, it is not possible to 
solve these equations using normal Gaussian quadrature, 
and special algorithms must be implemented. While 
these techniques have been outlined in the literature pre- 
viously p^ . Ilfil | , they have not been described in detail. 
In this appendix we will describe numerical routines used 
to calculate integrals of the form 



which results in the integral 



dfc F{\q-k\)G{k) 



(A.l) 



when F(k) and G(k) are only known on a grid of equally 
spaced wave vectors and and how to solve equations of 
the form 



Fq{t) = aF^(t) + / Au Mg{Fg, t - u)Fg{u) (A.2) 
Jo 

for many decades in time. The dot denotes a derivative 
with respect to time. It is important to note that the 
function Mq not only depends on time t, but also on the 
functions Fq{t). 

The integrals of the memory function in the mode- 
coupling theory are all of the form given by Eq. lA.ll To 
calculate the integral, the functions are generally only 
known on a grid of equally spaced wave vectors. Fur- 
thermore the calculation of these integrals are the most 
computationally expensive part of the program, and ex- 
trapolation between grid points can increase the calcula- 
tion time significantly. The first step is to introduce the 
change of variables p = k — q/2, then convert to spher- 
ical polar coordinates. This allows integration over one 
angular variable and the integral becomes 



27r / p sin (j) dp dcj) F 



G 



Then make another change of variables to 



y = 



= ^Jp"^ + + pqcosi 
~ \/'{P + 9^/4 ~ pq cos( 



(A.3) 

(A.4) 
(A.5) 



2n f°° f'^^'^ 

— / a;da; / ydyF{y)G(x) 

1 Jo J\x-q\ 



(A.6) 



which can easily be calculated numerically using any 
number of quadrature techniques. There is one technical 
issue in using Eq. IA.6I As q ^ the integral over y goes 
to zero, but the integral itself does not go to zero, see 
e.g. equation 1231 While the effect is negligible for larger 
wave vectors, for small q it is better to expand Eo. lA.ll in 
a Taylor series around q = 0. We found that this had to 
be done to obtain accurate small q values of varies quan- 
tities, e.g. the non-ergodicity parameter. In this work 
the Taylor series approximation was always used to cal- 
culate the integrals of the memory functions for the two 
smallest wave vectors. 

Most of the equations solved in this paper has the 
form given by Eg. I A. 21 and they must be solved for many 
decades in time. The basic algorithm is as follows. First 
an arbitrary time interval At is broken into AN equal 
segments of size 6t = At/4N. It is assumed that the 
value of Fq{t) for the first 2N segments is known. For 
each future time an equation of the form 



Fq{t)^aH{Fq,t)+b. 



(A.7) 



is solved. When Fq{ti) has been calculated for all AN 
times the time interval At' = 2/S.t is doubled. The new 
time integral in divided into 47V equal segments of size 
5t' = 25t. Since F{t) has been calculated up to Ai'/2, 
a mapping of {Fq{ti)} — > {Fq[tj)} can be defined where 
ti — iSt and tj ~ jSt' and thus we know Fq(t) for the first 
2A^ segments of the new time interval, and the procedure 
is repeated. 

First we will describe how to convert Eq. IA.2I into an 
equation of the form I A. 71 Break the integral in Eq. IA.2I 
into two integrals, then integrate the integral starting at 
t = by parts to get 



du Mq{t~u)Fq{u) = Mq{t-t2)Fq{t2)- Mq{t)Fq{0) 



duMqit - u)Fq{u) 



+ / duMq{t-u)Fq{u) (A.8 

Jt2 

Next make a change of variables in the second integral in 
Eg. lA. 81 to T = t — u and break both integrals in Eg. IA.8l 
into integrals of length At/2. This results in the following 
exact form of the integral in Eq. IA.8L 

Mq{t ~ t2)Fq{t2) 



Mq{t)Fq{0) 
"1 pt 



/ duMqit - u)Fq{u) 

/ ' duMq{u)Fq{t-u). 
.7 = 1 -^'j-l 



(A.9) 
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Use the approximation that 



1 



duA{u)B{u) w {A{tj) ~ A{tj^i)}j- I duB{u) 

« {A{t,)-A{t,^,)}I[B{t,)] (A.IO) 

where I[B{tj)] = (l/2)[B(tj_i) + B{tj)]. We approxi- 
mate Fq(ti) by 

FS) « - + (A.ll) 

but other approximations for the derivative can be used. 
Put all this together to get the equation 



(A.12) 



where 



C3 - —Fq{ti^l) - —Fq{ti^l 



G2 = X[F,(ti)]-F,(0) 
St 

+Mq{U^a)Fq{U2) - MqiU^i) 
-Fq{U^i)J [Mq{U2)] 

~ J2 {Mq{U-i) - Mq{U.,+^)}I[Fq{t,)] 

(A.13) 



which can be easily recast into the form of Ea. lA.7l Recall 
that Mq{ti) depends on Fq{ti). Equation IA.7I can be 
solved using any number of techniques used to find a 
fixed point of a set of equations. There is one equation 
of the form IA.7l for every q vector, and Mq{ti) depends 
on each Fq{ti). 

It remains to describe the mapping from At to At' — 
2At. For 1 < j < 2N, 



Fq{t2.,) 
Mq{t2.,) 



Mq{t,). 



(A. 14) 
(A.15) 



For 1 < J < A^, 



Q.h{I[Fq{t2^)]+I[Fq{t2^-l)]} ~^ I[Fqit,)] 

(A.16) 

{).h{l[Mq{t2,)]+l[Mq{t2^-l)]} ^ T[Mq{tj)] 

(A.17) 

For TV + 1 < j < 2iV 



-^{Fq{t2i) + ^Fq{t2^-l) + Fq{t2^-2)} ^ I[Fq{tj)] 

(A.18) 

\ {Mq{t2i) + AMq{t2^-l) + Mq{t2^-2)} ^ 2: [Mq{tj)] . 

(A.19) 

Now there are values for F, M, J[F], J[Af] from < t < 
At' /2 and they can be calculated for At' /2 < t < At' . 
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